/*
Plots of regression estimates for individual entry events
1 => Lipitor
9 => Crestor
10 => Vytorin
12 => Generic Zocor
*/


set more off
clear
use user_quarter_choice_cholesterol



// focus on cholesterol
tab indic_id
keep if indic_id == 5




// BASIC: LOOK AT IMPACT OF FIRST CHOICE ON LATER CHOICE 
egen first_q = min(quarter), by(enrolid)
gen first_choice_temp = generic_id if quarter == first_q
egen first_choice = max(first_choice_temp), by(enrolid)

gen first_generic_temp = generic if quarter == first_q
egen first_generic = max(first_generic_temp), by(enrolid)

gen first_xr_temp = xr if quarter == first_q
egen first_xr = max(first_xr_temp), by(enrolid)


gen period = quarter - first_q
gen chose_initial = (generic_id == first_choice) & (generic == first_generic) & (xr == first_xr)


drop first_choice_temp first_choice





/////////////////////////////////////////////
// IDENTIFIED ANALYSIS: USING ENTRY EVENTS //
/////////////////////////////////////////////



// need name for overlap analysis
merge m:1 generic_id using generic_name_mapping
keep if _merge == 3
drop _merge


rename generic_id curr_generic_id
rename gennme curr_gennme
rename generic curr_generic
rename xr curr_xr
drop drugid

// user is in the dataset only based on first use
joinby enrolid indic_id using "user_treatment_indicators_cholesterol"


// filter out irrelevant entries
// 1) not either treatment or control for that entry event
keep if treat == 1 | control == 1
replace treat = 0 if control == 1

//replace treat_alt = 0 if control_alt == 1

// 2) quarter came before entry quarter
gen entry_quarter = qofd(entry_date)
gen rel_quart = quarter - entry_quarter
//keep if rel_quart >= 0

// A) generate the outcomes
gen chose_entry_drug = (curr_generic_id == generic_id) & (curr_generic == generic) & (curr_xr == xr)
tab chose_entry_drug


// B) analyze
gen treat_generic = treat * generic 
gen treat_xr = treat * xr
// is stickiness bigger for generic choice? (doesn't have to be)


// Does this need correcting?
gen initial = chose_entry_drug if quarter == first_q  // rel_quart == 0 (old code; this is probably cleaner; use their first choice, otherwise some not there in period 0)
egen initial_choice = max(initial), by(enrolid entry_num) // unnecessary? maybe only 1 entry_num per enrolid


gen initial_generic_temp = generic if quarter == first_q
egen initial_generic_user = max(initial_generic_temp), by(enrolid)


/////////////////////
// plot ols vs. IV //
/////////////////////
// http://www.stata.com/meeting/germany14/abstracts/materials/de14_jann.pdf, slide 27


// NEW: add coefplot
matrix ols = J(18, 3, .)
matrix coln ols = coeff l95 u95
matrix rown ols = 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

rename initial_choice initial_choice_temp
gen initial_choice = 1

forvalues i=3(1)20 {
	display "Offset: `i' quarters"
	reg chose_initial initial_choice if period == `i' & initial_generic_user == 0, nocon r
	matrix b= e(b)
	matrix v = e(V)
	
	matrix ols[`i'-2,1] = b[1,1], b[1,1] - sqrt(v[1,1]) * 1.96, b[1,1] + sqrt(v[1,1]) * 1.96
	
}

replace initial_choice = initial_choice_temp

matrix iv = J(18, 3, .)
matrix coln iv = coeff l95 u95
matrix rown iv = 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

forvalues i=3(1)20 {
	display "Offset: `i' quarters"
	ivregress 2sls chose_entry_drug (initial_choice=treat) if rel_quart == `i' & generic == 0, r

	matrix b= e(b)
	matrix v = e(V)
	
	matrix iv[`i'-2,1] = b[1,1], b[1,1] - sqrt(v[1,1]) * 1.96, b[1,1] + sqrt(v[1,1]) * 1.96
}




/*
Do it by entry event
*/

foreach g in 1 9 10 12 {

preserve

keep if entry_num == `g' // look at branded entries


matrix iv = J(18, 3, .)
matrix coln iv = coeff l95 u95
matrix rown iv = 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

forvalues i=3(1)20 {
	display "Offset: `i' quarters"
	ivregress 2sls chose_entry_drug (initial_choice=treat) if rel_quart == `i', r

	matrix b= e(b)
	matrix v = e(V)
	
	matrix iv[`i'-2,1] = b[1,1], b[1,1] - sqrt(v[1,1]) * 1.96, b[1,1] + sqrt(v[1,1]) * 1.96
}

coefplot (matrix(iv[,1]), ci((iv[,2] iv[,3]))), ///
vertical nooff xtitle("Quarter") ytitle("Estimate (Fraction of Users)") yscale(range(0 1)) ylabel(0 0.25 0.5 0.75 1, grid)


graph save "coefplot_iv_entry_`g'", replace
graph export "coefplot_iv_entry_`g'.pdf", replace

restore

}